mtmM.. 


m^2^-wm 


Computer  Science  Department 


TECHNICAL  REPORT 


An  Additive  Schwarz  Method  for  the  p-version 
Finite  Element  Method 


Luca  F.  Pavarino  * 


Technical  Report  580 

September  1991 


■  «*»ui^ 

'~^      bo; 

o      -C  •*-• 

'-'        4J   -H 

0)  c 

E  -H 

«W 

N 

O          U    c 

CO         rn   o 

LD   t,    5  -H      . 

NEW  YORK  UNIVERSITY 


I 

El    O 


U 

a.' 

o 
o 


U    S-)    O 

>  -P 

I    <u 
a,  c 


> 


4-1    0)  4-) 

•H  ^  C 

t3  -P  Q) 

(0    k-j  0) 
O 


Department  of  Computer  Science 
Courant  Institute  of  Mathematical  Sciences 

251  MERCER  STREET,  NEW  YORK,  N.Y.  10012 


An  Additive  Schwarz  Method  for  the  p-version 
Finite  Element  Method 


Luca  F.  Pavanno  * 
Technical  Report  580 

September  1991 


*  Courant  Institute  of  Mathematical  Sciences,  251  Mercer  Street,  New  York,  NY 
10012.  This  work  was  supported  in  part  by  the  Applied  Math.  Sci.  Program  of  the 
U.S.  Department  of  Energy  under  contract  DE-FG02-88ER25053  and,  in  part,  by  the 
National  Science  Foundation  under  Grant  NSF-CCR-8903003. 


Abstract 

The  additive  Schwarz  method  was  originally  proposed  for  the  h- 
version  finite  element  method  for  elliptic  problems.  In  this  paper,  we 
apply  it  to  the  p-version,  in  which  increased  accuracy  is  achieved  by 
increasing  the  degree  of  the  elements  while  the  mesh  is  fixed.  We  obtain 
a  constant  bound,  independent  of  p,  for  the  condition  number  of  the 
iteration  operator  in  two  and  three  dimensions.  The  result  holds  for 
linear,  self  adjoint,  second  order  elliptic  problems  and  for  quadrilateral 
elements. 


1  Introduction. 

In  the  p-version  of  the  finite  element  method,  the  degree  of  the  piecewise 
polynomial  elements  is  increased  in  order  to  achieve  the  desired  accuracy, 
while  the  mesh  is  fixed.  This  is  in  contrast  to  the  standard  h-version  where 
fixed  low  order  polynomial  elements  are  used  and  the  mesh  is  refined  in  order 
to  obtain  accuracy.  For  an  overview  and  basic  results  about  the  p-version, 
see  Babuska  and  Suri  [1].  In  this  paper,  we  study  a  domain  decomposition 
method  using  p-version  finite  elements  in  the  framework  provided  by  the  ad- 
ditive Schwarz  method  (ASM).  We  consider  linear,  self-adjoint,  second  order 
elliptic  problems  and  brick-shaped  elements  in  the  finite  element  discretiza- 
tion. We  show  that  the  condition  number  of  the  ASM  iteration  operator  is 
bounded  by  a  constant  independent  of  p.  The  proof  is  similar  to  the  one  in 
Dryja  and  Widlund  [4]  and  is  based  on  Lions'  partitioning  lemma. 

This  paper  is  organized  as  follows.  In  Section  2,  we  define  a  model 
problem  and  introduce  its  discretization  with  the  p-version  finite  element 
method.  In  Section  3,  we  review  the  basic  framework  of  the  ASM  and  apply 
it  to  our  model  problem  with  square  elements  in  two  dimensions.  An  ASM 
using  brick-shaped  elements  in  three  dimensions  is  considered  in  Section  4. 
For  an  example  of  a  method  similar  to  ours,  but  using  the  h-version  finite 
element  method,  see  Bramble  et  al.,[2]. 

2  The  Model  Problem. 

We  consider  a  model  problem  for  linear,  self  adjoint,  second  order  elliptic 
problems,  on  a  bounded  Lipschitz  region  Q.  The  discrete  problem  is  given 
by  the  p-version  finite  element  method.  For  simplicity,  we  first  consider  the 
following  problem  in  R^: 

-Au     =     f  in  fi, 

u    =    0  on  dCt . 

The  standard  variational  formulation  of  this  problem  is  : 
Find  u  e  V  =  H^iQ)  such  that 

a(u,v)  =  f{v),  V  uG    F, 

where  the  bilinear  form 


a(u,v)  =   /  Vu  •  Vv  dx 
Jn 


defines  a  semi-norm  |w|//i(fi)  =  (a{u,u))^/'^  in  H\Q),  and  a  norm  in 
V  =  H^{il).  Our  analysis  works  equally  well  for  any  general  self  adjoint, 
continuous,  coercive,  bilinear  form  6(-,-)»  since  the  H^  norm  and  the  one 
induced  by  6(-,  •)  are  equivalent: 

c|u||^j  <  b{u,u)  <  C|u|^i  . 

A  triangulation  of  the  region  Q.  is  introduced  by  dividing  it  into  non- 
overlapping  brick-like  elements  fl,  ,i  =  l,---,Ne  .  We  suppose  that  the 
original  region  is  a  union  of  such  elements  and  we  denote  the  mesh  size  by 
H. 

We  define  Qp  to  be  the  set  of  polynomials  of  degree  less  then  or  equal 
to  p  in  each  variable,  i.e.  in  two  dimensions 

Qp  =  span{x'y^  :  0  <i,j  <p} 

and  we  discretize  the  problem  with  continuous,  piecewise,  degree  p  polyno- 
micil  finite  elements: 

yp  =  {</,  e  c°{il) :  (Alfi,  G  Qp,  i  =  1,  •  •  • ,  N^} . 

Then  the  discrete  problem  takes  the  form  : 
Find  Up  e  VP  such  that 

a{up,Vp)  =  f{vp),  V  Vp  G    v.  (1) 

For  simplicity,  we  will  analyze  square  and  brick-shaped  elements,  but,  using 
affine  mappings  onto  the  reference  square  and  cube,  our  analysis  works  also 
for  general  quadrilateral  elements. 

3      An  additive  Schwarz  method  using  square  ele- 
ments in  two  dimensions. 

The  additive  Schwarz  method  was  originally  developed  for  the  standard  h- 
version  finite  element  method  and  we  refer  the  reader  to  Dryja  and  Widlund 
[4]  or  Dryja  [3]  for  more  detailed  discussions.  We  work  now  in  dimension 
two  and  with  square  elements  fi,  . 

Let  N  be  the  number  of  interior  nodes.     Our  finite  element  space  is 
represented  as  the  sum  of  N-l-1  subspaces 


The  first  space  Vq   serves  the  same  purpose  as  the  coarse  space  in  the  h- 

version.  Here: 

Vq"  =  V^,  i.e.  the  space  of  continuous,  piecewise  Qi  functions  on  the  mesh 

defined  by  the  elements  fi,; 

V.P  =  VPn^o(^i)   where  fij  is  the  2H  X  2H  open  sqiiare  centered_  at  the 

i-th  vertex.    In  other  words  il'-  is  the  interior  of  fi,,  [Jili-^  U^is  U^m  »  see 

figure  1  below.  As  in  the  h-version,  the  algorithm  consists  in  solving,  by  an 


fi„ 


i^M 


n„ 


i)„ 


Figure  1:  the  substructure  fi, 
iterative  method,  the  equation 

Pup  =  (Po  +  Pi  +  ...  +  Pn)Up  =  g'p  , 
where  the  projections  P,  :  V^  —*  Vf  are  defined  by 

a{P,Vp,4>p)  =  a{vp,(t>p)  ,  'i(f>peVf. 


(2) 


(3) 


The  following  is  the  main  result  of  this  paper.   The  proof  is  first  given 
for  two  dimensions  and  then  extended  to  three  in  Section  4. 

Theorem  1  The  operator  P  of  the  additive  algorithm  defined  by  the  spaces 
V-^  satisfies  the  estimate  k(P)  <  const,  independent  of  p. 

Proof.  The  proof  is  similar  to  a  result  given  in  Dryja  and  Widlund  [4]  for 
the  h-version.  A  constant  upper  bound  for  the  spectrum  of  P  is  obtained 
by  noting  that  for  i  >  1 


a{PiUp,Up)  =  a{P,Up,P,Up)  -  a^:{P,Up,  P,Up)  <  aQ,^[Up,Up). 


Each  point  is  covered  by  no  more  than  four  subregions  Q',  and  the  norm  of 
Po  is  equal  to  one,  therefore  Xmax  <  5  . 

A  lower  bound  is  obtained  by  using  Lions'  lemma,  see  Lions  [6]  for  the 
case  A^  =  2;  a  proof  is  also  given  in  Widlund  [8]. 

Lemma  1  Let  Up  -  Eilo^P.^  '^^^^  "?,•  e  K,  ^  a  representation  of  an 
element  o/ V  =  Vb  +  Vi  + \-  Vn-  U  the  representation  can  be  chosen  so 

that 

N 

Y^  a{Up^„  Up,,)  <  Coa{up,  Up)  ,  Vup  G  V  , 

1=0 

then 

^min(  ')  ^  ^0 

We  have  to  define  a  partition  of  the  finite  element  function  Up  and  obtain 
a  good  bound  of  Cq.  We  know  from  Strang  [7]  that  there  exists  a  linear  map 
ii   :  yp  — >  V',  which  satisfies 

ll"p  -  Awp||i2(n)  <  ClH^\up\]f^^^^  (4) 

and 

\Up    -    /l«pi//l(Q)     <    C2\Up\fjlf^Qy  (5) 

Let 

Up,o  =  hup,  Wp  =  Up  -  Upfi  . 

In  order  to  define  Up,  G  Vj' ,  we  consider  a  particular  partition  of  unity  {^,} 
consisting  of  the  standard  basis  functions  for  Qi{=  Vi)  : 

N 

OieVi,        supp(e,)  ^Q[,        o<o,<i,        ^  0,(x,  y)  =  i . 

i=l 

Now  0iWp  is  an  element  of  Qp+i  vanishing  outside  Q[.  Since  in  our  partition 
we  need  an  element  of  V,''  ,  we  interpolate  back  0,Wp  into  V^^  .  We  define 
this  interpolation  operator  Ip  on  one  of  the  four  elements  il,  of  il[]  on  the 
other  three  it  is  completely  analogous.  We  transform  this  element  into  the 
reference  square  [—1, 1]  X  [-1, 1]  .  S,  is  1  at  one  vertex  of  fi,  and  0  at  the 
other  three,  so  it  has  one  of  the  four  forms 

e,=  ^-{x±l)iy±l). 


We  define  Up,,  =  Ip{0tWp)  as  the  polynomial  in  Qp  interpolating  9iWp  at  the 
(p  +  1)^   points  (Xn»^m)i  where  the  x'^s  are  the  zeros  of  the  polynomial 

Cp+i{x)  =   /     Lp(s)ds  .  (6) 

Here  Lp(s)  is  the  Legendre  polynomial  of  degree  p.  This  definition  makes 
sense  for  p  >  1  because  £p+i  has  p+1  distinct  real  zeros  in  [—1,1].  In  fact 
£p+i(±l)  =  0  and  p—  1  roots  interleave  those  of  Lp,  which,  as  is  well  known, 
has  p  distinct  real  zeros  in  [-1, 1].  We  define  Cq  —  1.  We  remark  that  even 
if  this  definition  of  the  interpolation  operator  is  local,  we  obtain  an  element 
in  Vf  =  V^  f]  //d(nj) .  In  fact,  Up^i  is  continuous  across  element  boundaries, 
because  on  the  boundary  of  each  element  we  have  p+1  interpolation  points 
and  Up,  is  a  polynomial  of  degree  p.  Since  Ip  is  a  linear  operator,  we  have 

N 

X]  "p.'  =  "p  -  "p,0  • 
:  =  1 

We  note  that 

/plg^  =  identity  . 

Since  |  •  \ffi    is  a  semi-norm  on  Qp  ,  it  is  natural  to  introduce  the  quotient 
space  Qp  =  Qp/Qo  ,  on  which  |  •  |//i    is  a  norm.  Clearly  dimQp  =  (p  +  1)^, 
while  dimQp  =  (p  +  1)^  -  1  . 
We  now  establish: 

Lemma  2  The  interpolation  operator  Ip  :  Qp+i([-l,  1]^)  ^  Qp([-1, 1]^)  is 
uniformly  bounded  in  the  \  ■  |//i    norm,  i.e. 

|/p(/)l//i  <  const.\f\„^  ,  V/  G  Qp+i([-l,  I?)  . 

Proof.  If  /  is  a  function  of  x  only,  then  Ipf  is  a  function  of  x  only  and 

\ipif)\h  _  Wj-MDWh 

l/"P  11^112 

Similarly,  if  /  is  a  function  of  y  only 

\ipim.  _  iit^p(/)iii3 

!/!?/■  "  m\h  ■ 


In  general,  if  both  1^  and  |^  do  not  vanish,  it  is  easy  to  see  that 


\W)\h     _    \\fMf)\\h  +  \\fyW)\\h 


< 


-  \\91\\2  "^  \\dl\\2  ■  *-    -* 

We  consider  the  two  terms  separately.  By  symmetry,  we  need  only  to  study 
the  first  term.  We  form  a  basis  for  Qp+i  using  the  polynomials 


1     C,(y) 

v^II-c.IIl^' 


</>o.j(x,y)=     ^ii/ii       ,  1<;<P+1, 


and 


1       C,(x) 
<l>i,o{x,y)=-7=../\       ,  l<i<P+l, 

<i>.,(^,y)=..^'^"j    ..y.f^  ,  l<t,j<p+l.  (8) 

We  can  disregard  the  space  spanned  by  the  <Ao,j's  because  every  function  /o 
in  that  space  will  not  contribute  to  the  x-term  we  are  considering.  In  fact, 
if  /  =  /o  +  /,  we  have 

In  the  resulting  space  Qp+i  of  dimension  p'  =  (p  +  2)(p  +  1),  we  choose  the 
order 

<f>\Oi<f>ll,  ■  •  ■  ,'t>\p+l,4>20,'l>2\,  ■  •  ■  ,4>2p+l,  •  •  •  T  <f>p+l  ,0t  4>p+l  ,1 ,  ■  ■  ■  ,<t>p+l,p+l- 

If,  for  simplicity,  we  relabel  the  basis  as  {(l>k{x,  y) ,  k  =  1,  2,  ■  •  •  ,p'}  ,  we  have 


f{x,y)=  ^afc<?iit(x,y), 
k=i 

rrl  I  <Pk    II  <t>k  &Qp  ■ 


The  last  relations  follow  from  the  fact  that  the  coordinates  of  the  inter- 
polation nodes  (a;„,a;m)  are  the  zeros  of  £p+i  and  /pig  is  the  identity. 
Hence 


i^/m2        /v^      9(t>k  ^     d<l)i.         ,--         ,d<t>k   d(t>i 


dx 


k=i    ^^';^  '^^ 


k,l 


—  a   SxOt , 

where  (5x)*:,/  =  (^^?  ^)l'^  '^  ^  symmetric,  positive  definite  matrix  of  order 
p' .  Similarly, 


where 


j-^IvU)\\h  =  o^'^BS^Ba  , 


/  / 


B  = 


\ 


0 


is  a  projection  matrix  onto  Qp  .  Here  the  identity  matrices  are  all  of  order 
p  +  1,  the  last  zero  block  is  of  order  p  +  2  and  the  other  zeros  are  of  order 
1.  Therefore, 

WfMDWh       a^'BSxBa 

The  proof  of  the  lemma  follows  from  a  bound  for  the  eigenvalues  A  of  the 
generalized  eigenvalue  problem 


BSxBa  =  XSxO  . 


(9) 


The  structure  of  S^  can  be  made  explicit  by  using  the  orthogonality  of  the 
Legendre  polynomials  and  the  formula 


I\  ^"^'^ 


ds 


2n  +  1 


n  >  1 


(10) 


In  fact,  if  M^,y)  =  y^^   and  </>,(x,y)  =  ^fffj  ,  we  find 


{s.hi  =  i^M^.y)-J'(-^y))Li,  =  (||L._i||  ||£,||'|iIZj  \\Cm\\  ^^- 


d_ 
dx 


d_ 
dx 


L,.i{x)   Ln-i{x)        Cjjy)   Cmjy). 

This  expression  differs  from  zero  iff 

'j-2 
n  =  i  and  m  =  <       j        , 

and  therefore  each  row  of  S^  has  at  most  three  nonzero  elements.  They  are 


(5,)fc/  =  < 


c,      c 


I      '->  ''T  +  2      \ 

(    MKj||'|K;+2ir'^J 

The  only  exceptions  to  this  rule  occur  when  one  of  the  indexes  is  0, 1  or  2. 
In  this  case  we  cannot  use  formula  (10),  but  we  can  use  £o  =  1  =  ^o  and 
Ci{x)  =  /f  J  (is  =  X  -\-  I  -  Li  +  Lq.  The  exceptional  elements  are  then 

d  d  1     jCi(y)^  /3 

(5.),,  =  (^</'.,o,^</'.,i)  =  1  •  ^Tf '  1^)^?  =  V  2  =  '°  • 


(5.),,  =  (-<^,o,^</>.,2)=l-(^,^|^ 
r?^        r  ^A      ^^    ^      1    r-^i'^/)  ^2(2/).  /5 


^j  =  -V6  =  ^°- 


Cl 


This  shows  that  Sx  has  the  structure 


/  A 


p+2 


h+2 


h+2  ) 


Each  block  Ap+2  is  symmetric,  pentadiagonal  and  of  order  p  +  2 


Ap+2   = 


/I       Co      60 

Co        1        Ci  61 

60     Ci      1  0     62 

61      0  1      0 

62  0      1 


\ 


\ 


bp-2 

1  0       6p_i 


6p-2       0  1 


'Jp-2 


0 


bp-1        0  1      / 


The  elements  co,ci  and  bo,bi  have  been  defined  above  and 

By  using  (10)  and  ||Z'n|||2r_i  ji  -  2^^  ,  we  can  compute  the  6j's  explicitly 

4 


(£_,, £j+2)  -  - 
and  therefore 


(2j-3)(2;-l)(2j  +  l)' 

2 

~(2j-l)(2j-fl)(2j  +  3)' 


j>2, 


IA||  =  ^/3, 


bj 


1      (2j-3)(2;  +  5) 


j>l,  (A,  A) 


J>2 


2V(2i-l)(2;  +  3)' 
Our  generalized  eigenvalue  problem  (9)  can  now  be  written  as 

/   Bp+2 


(11) 


B 


p+2 


"2 
"3 


B 


P+2 


0; 


/  A 


=  X 


P+2 


lp+2 


V  "p+1  / 

\  /  «1  \ 

"2 
03 


lp+2 


^P+2    ) 


\    «P+1    / 


where 

B„4.2  = 


P+2   -    I  0 


The  last  zero  is  a  scalar  and  the  q^'s  are  vectors  of  length  p  +  2.    This  is 
equivalent  to  Op+i  =  0  and 

I    ^"^^     n    I  ".  =  AAp+2a.-  ,  l<i<p.  (12) 


0 


But 


Ap+2  -    [       i,T         1 


with  6^  =  (0,  •  •  • ,  0, 6p_i,  0)  =  bp_\e^  ,  where  e  is  the  p  —  th  column  of  the 
identity  matrix  of  order  p  +  1.  Therefore,  (12)  is  equivalent  to 

^p+i  ^\  ( "^  \  - x( ^p+i  M  r^' 

0      0  )[w  )  -^[    b-^     1/V^ 

I.e. 

Ap+iv  =  X{Ap+i-bl_iee'^)v.  (13) 

Since  ee-^  has  rank  1,  we  see  immediately  that  we  have  p  eigenvalues  equal 
to  1,  corresponding  to  the  eigenvectors  v  with  Vp  =  0  .In  order  to  find  the 
only  non-trivial  eigenvalue,  we  apply  the  Shermann- Morrison  formula  (see 
Golub  and  Van  Loan,  [5],  pg.  51) 

{A  +  uv'^)-^  =  A-^ ^r-—-A-^uv'^A-^  (14) 

l  +  v^A-^u  ^      ' 

to  >lp+i  —  6p_iee^  and  obtain 


J\-\  _  4-1     I  "p-i  1-1  ^„r^-i 


(^p+i  -  ^p-i^^'')-'  =  Kl^  +  7-T2-T5^7^r-^p-+i^e^  Vi 


With  this  formula,  we  have  reduced  (13)  to  the  standard  eigenvalue  problem 
Mv  =  (Ap+i  -  b\_^et^ Y"^ Ap^xv  -  Xv  ,  (15) 

with 

M  -  I  A 2zi 4-1   ppT 

^        °p-l^    ^p+1^ 

10 


Again  we  can  see  that  the  eigenspace  corresponding  to  A  =  1  has  dimension 
p  by  substituting  v  Gnullspace(^p^ie)  into  (15).  We  obtain  the  non-trivial 
eigenvalue  by  substituting  v  =  -^p+i^  iiito  (15): 

with 

A^i+/y^/;!:  =     ^  ,  ■      (16) 

In  order  to  obtain  an  upper  bound  for  A,  we  need  some  properties  of  A~_^_i 
and  the  sequence  of  6j's  .  Since  from  (11),  6p_i  <  \,  we  find  that 


A< 


1  _  IpT 4-1   p 


What  remains  is  to  find  a  bound  on  e^ Ap^^e  =  the  p— th  diagonal  element 
of  ^p+i  •  Let  Cp  =  det{Ap) .  By  Cramer's  rule, 

e^A;l,e=^,  (17) 

flp+i 

where  the  cofactor  dp  is  the  determinant  of  the  matrix  obtained  from  Ap+i 
by  deleting  the  p-th  row  and  column.  Applying  Laplace's  theorem  for  the 
expansion  of  determinants  to  Op+i  and  dp  ,  it  is  easy  to  prove  the  following 
recurrence  relations: 

flp+i  =  flp  -  6p_2(ap_2  -  6^_3ap_3) ,  p  >  4,  (18) 

flp  =  flp-l  -  6p_2«p-2  ,  P>4.  (19) 

(18)  can  be  written  as 

ap-2  -  6p_3ap_3  ,  (20) 


iff  ap_2  >  6p_3ap-3  •  (21) 


ap- 

flp+i 
-2 

which 

shows 

that 

ap 

> 

flp+i 

11 


Now,  Ap+i  is  positive  definite,  because  a^ Ap^ia  defines  a  L^  norm  of  a 
function  with  components  {a,}  .  Therefore  A~^i  is  positive  definite.  Sub- 
stituting (19)  into  (17)  and  using  (20),  we  get 

0  <  eA-l,e  =  °^-^  "  ^'-^""-^  =  ";r  '  "'^"  =  7J-(1  -  ^^)  •      (22) 

'^^  flp+l  Op-lOp+l  Op-l  «P+1 

This  implies  that  1  -  ^^  >  0  ,  i.e.  Cp+j  >  ap+2  ,Vp.  By  (21)  we  then  have 
ap+2  >  ^p+iflp+i  ,  i-e-  f^  >  ^p+i  ,Vp.  Hence 

e%"+i^<IT-(l-''p+i)- 
"p-i 

Since  limp^oo  ''p  =  4  S'Hd  6^  <  -j  ,  for  every  e  >  0  we  have  6p  >  ^(1  -  e)  for 
p  large  enough.  This  gives  us  the  desired  bound 

e%-ii'^<l^(l-^(l-0)  =  f^  =  3  +  .' 

for  p  large  enough,  and  finally 

1  A 

A  < 


I-i(3  +  €')       1-e' 

In  other  words,  A  <  const,   uniformly  in  p . 

Numerical  experiments  in  MATLAB  show  that  a  stronger  result  is  actually 

true  :  limp_oo  A  =  2  . 

In  conclusion,  we  have  found  a  bound  for  the  x-term  in  (7): 

\\£ip(f)\\h  .      , 

sup         „, '^  <  const.  . 


Reasoning  in  the  same  way  for  the  y-term,  we  find  from  formula  (7) 


— ,  .,0  "     <  const. 


We  can  now  conclude  the  proof  of  the  theorem  by  applying  Lemma  2  to 
bound  the  //'  norm  of  each  component  Up^,  =  Ip(6iWp)  over  a  single  element 
ilk. 
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On  a  square  element  Qk   of  side  H  ,  \^\   and  |^|   are  bounded  hy  l/H 
and,  by  construction,  ||^t||L°°  <  1  •  Therefore 

2 

Since  at  most  4  components  Up,,  are  nonzero  for  any  element  ilk,  we  obtain, 
when  summing  over  i, 

^2  2  2 

1=1 

and  summing  over  all  the  elements  Q.k 

^2  2  2 


Using  equations  (4)  and  (5),  we  can  conclude 

i(Q)  <  8C(2Ci  +  C2)|itp|^,(j2j  =  cuusi  |"p|^i(Qj 


Zl  l"p--l//>(fi)  -  8'^(2Ci  +  C'2)|«p|^wj))  =  const  lupl^i 


1=1 


4     An  Additive  Schwarz  Method  using  brick-shaped 
elements  in  three  dimensions. 

Theorem  1  can  be  extended  to  dimension  three  and  brick-shaped  elements 
and,  by  induction,  to  an  arbitrary  dimension.  We  define  in  this  case 

Qp  =  span{x' y^ z^   :  0  <  i,j,k  <  p} 

and  we  assume  that  the  region  fi  C  -ft"'  is  the  union  of  non-overlapping 
brick-shaped  elements  fi,  of  side  H .  If  A'^  is  the  number  of  interior  nodes, 
we  represent  V^  as 

l/P  =  Fo"  +  V^i"  +  •  •  •  +  V^^  , 
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where  V^  =  V^  and  V;"  =  V^  n  H^{n',)  .  Now  Q;  is  the  open  cube  of  side 
2H  centered  at  the  i-th  interior  node.  We  now  prove  the  main  theorem  in 
dimension  three. 

Using  the  same  notation  as  in  the  two  dimensional  case,  we  use  Lions'  lemma 
to  bound  the  spectrum  from  below.  The  upper  bound  is  obtained  as  before. 
The  partition  of  the  finite  element  function  Up,  required  by  Lions'  lemma,  is 
constructed  using  again  the  results  of  Strang  [7]  to  find  Up,o  and  a  piecewise 
linear  partition  of  unity  {^,}.  On  the  reference  cube  [—1,1]'',  ^,  can  have 
one  of  the  eight  forms 

0,  =  ^{x±l)iy±l){z±l). 

Recalling  that  Wp  =  Up  -  Up^  ,  we  define  Wp,i  =  Ip{9iWp)  as  the  polynomial 
in  Qp  interpolating  OiWp  at  the  {p  +  1)^  points  {xi,Xm,Xn)  ,  where  the  a;„'s 
are  the  zeros  of  the  integrals  of  the  Legendre  polynomials  £p+i  defined  in 
(6).  Working  again  with  the  quotient  spaces  Qp  —  Qp/Qo  ,  we  can  establish 

Lemma  3  The  interpolation  operator  Ip  :  Qp+i([-l,  1]'^)  — *  <3p([— 1, 1]"')  is 
-imformly  hounded  n;  the  |  •  |^i  norm,  i.e. 

'Mf)\H^  <  const.\f\H,  ,  V/  G  gp+i([-l,  If)  . 

Proof.  Since 

\ipif)\l^  ^  Wfjpimh  ,  lll;^p(/)lli^  ,  \\i-Mf)\\l 

\f\l,        -         \\dl\\2  "^         ||9/||2  +         |iaX||2  '  ^^'^f 

i-'i//'  llaxllL2  II  ay  I1l2  llazllL2 

we  need  only  to  consider  one  of  the  three  terms,  for  example  the  x-term.  As 
before,  if  in  one  term  the  denominator  is  zero,  we  can  prove  that  (23)  is  still 
valid  after  dropping  that  term.  We  form  a  basis  using  the  polynomials 

and  the  same  modification  with  the  constant  term  4t  when  one  index  is 
zero;  cf.{8).  We  do  not  consider  polynomials  that  do  not  depend  on  x, 
because  they  do  not  contribute  to  the  x-term  we  are  considering.  We  now 
want  to  order  these  remaining  basis  functions  in  such  a  way  that  we  can 
use  the  results  obtained  in  the  two  dimensional  case.  For  every  fixed  /,  we 
have  a  two  dimensional  subspace.  We  order  the  basis  of  this  subspace  as  in 
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the  two  dimensional  case.  We  do  this  for  /  =  0, 1,  •  •  ■  ,P  +  1  ■  Relabeling  the 
^ij_i's  as  a  one  dimensional  array,  we  have 


I-— 1|2    ^  a^CxOt, 


d 


^W)\\h  =  (^' DCxDa 


We  are  interested  in  an  upper  bound  for  the  eigenvalues  of  the  generalized 
eigenvalue  problem 

DCrDa  =  Xda  .  (24) 

Here  the  stiffness  matrix  d,  of  order  p'  =  (p  +  2)^(p  +  1),  has  the  structure 


\ 


Ot-     — 


/    Sx 

cqSx 

bo 

Co  Si 

Si 

ciSi 

biSi 

boSx 

ciSi 

Si 

0 

b2Si 

biSi 

0 

Si 

0 

b2Si 

0 

Si 

\ 


6p_i5x      0 


bp-iSi 
0 


and  the  interpolation  matrix  Z),  of  the  same  size,  has  the  structure 


/  B 


D 


\ 


B 


B 


B 


\ 


0/ 


The  matrices  Si  and  B  have  been  defined  in  Section  3.   Clearly,  DCxD  is 
the  matrix 


'       Si        Co  Si      bo  Si 

CoSi        Si        CiSx  biSx 

boSx      CiSx         Si  0             b2Sx 

biSi       0  Si           0 

02  Sx  0                Si 


\ 


bp_iSx 

0 

bp_xSi      0  Si      J 
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where  S^  =  BS^B  .  Therefore  equation  (24)  gives  us  the  block  equations 
Sxioti  +  coa2  +  boas)  =.  XSj:{ai  +  coa2  +  boOs) , 

Sxibi-20it-2  +  oii  +  6,a,+2)  =  XSx{bi_20ti_2  +  a,  +  6,a,+2)  , 

5x(6p_iap_i  +  ckp+i)  =  A5x(6p_iap_i  +  Op+i) . 
These  are  all  generalized  eigenvalue  problems  of  the  form 

SrV  =  XSxV  ,  (25) 

where  u  is  a  linear  combination  of  some  a, .  But  this  is  the  same  generalized 
eigenvalue  problem  considered  in  the  two  dimensional  case;  see  eq.  (9).  We 
can  then  apply  the  two  dimensional  result  and  conclude  that  the  eigenvalues 
A  are  bounded  by  a  constant  independent  of  p.  Reasoning  in  the  same  way 
for  the  terms  in  y  and  z,  we  complete  the  proof  of  lemma  3. 

D 

In  order  to  complete  the  proof  of  theorem  2,  we  just  repeat  some  of  the 
arguments  given  in  the  two  dimensional  case.  We  first  apply  lemma  3  to 
bound  the  H^  norm  of  each  component  Up,,  =  Ip{9i'Wp)  over  a  single  element 
ilk-  We  then  sum  over  i,  recalling  that  at  most  8  components  Up^i  are  nonzero 
for  any  element  Qk  and  finally  we  sum  over  all  elements.  We  conclude  the 
proof  using  equations  (4)  and  (5). 

D 

Remark.  The  result  can  be  extended  to  any  dimension  by  induction.  The 
only  nontrivial  part  is  the  proof  of  the  lemma  about  the  interpolation  oper- 
ator. The  induction  step  from  dimension  n  to  n  +  1  is  essentially  analogous 
to  the  arguments  in  the  proof  of  Lemma  3.  We  consider  one  term  at  a  time 
and  order  the  basis  in  the  following  way:  first  fix  the  {n  +  l)-th  index  to 
be  equal  to  1  and  order  the  resulting  subspace  in  the  same  way  as  in  the 
case  of  n  variables;  then  fix  the  {n+  l)-th  index  to  be  equal  to  2  and  repeat 
the  process,  until  the  (n  +  l)-th  index  is  equal  to  p  +  1.  With  this  choice, 
the  stiffness  and  projection  matrices  have  a  block  structure  that  allow  us 
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to  reduce  the  (n  +  1)  dimensional  generalized  eigenvalue  problem  to  one  of 
dimension  n. 
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